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Abstract 

The combination of Finite Element Method (FEM) simulation and experimental photo- 
elasticity provides both qualitative and quantitative information about the stress field in a 
polymer composite and particularly along the fibre-matrix boundary. Investigations were 
made using model specimens containing up to five parallel glass fibres loaded at angles 
from 0° to 15°. The material properties of the geometrically equivalent FE models are 
calibrated using the experimental photoelastic outputs. In contrast to the experimental re- 
sults, FE simulations provide full 3D stress fields. For verification purposes, the 3D stress 
fields are reduced to two-dimensional synthetic photoelastic phase images, showing good 
agreement. Furthermore, detailed studies on the components of the stress tensor, partic- 
ularly statements concerning the shear transmitted by fibre-matrix boundary, progressive 
fibre failure and the effect of load angle variation are presented. 

Key words: A. Model material; B. Damage propagation; C. Finite element analysis 
(FEA); D. Photoelasticity, LSM. 



1 Introduction 



The progressive failure of long fibre composites, loaded in or at moderate angles to 
fibre direction, goes along with an increasing number of fibre fractures. The result- 
ing load redistribution between fibre and matrix depends on the fibre and matrix 
properties and the behaviour of the fibre-matrix interface. The interface is mainly 
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characterized by the transmitted shear which is itself strongly influenced by the 
load orientation [1,2,3]. Assumptions for correct load redistribution from broken to 
unbroken fibres are of special interest for the realistic description of materials us- 
ing models like fibre bundle models [1,4]. Especially the form of the load sharing 
relation determines the functional behaviour of the model in the breakup process. 

The breakup process can be observed experimentally with the non-destructive method 
of photoelasticity, using a confocal laser scanning microscope (LSM) [5,6]. Com- 
plex stress states are visualized by transillumination leading to 2D phase images 
clearly showing lines of equal difference between principal stresses (isochromates) 
and their gradients [7]. Unlike standard microscopes, the LSM provides informa- 
tion within any particular volume element of a transparent specimen while other 
information is effectively removed by a spatial filter [8]. Unfortunately phase im- 
ages are difficult to interpret, since in situ observation neither allows investigations 
on the stress variation in direction of specimen thickness nor on the composition 
of the stress tensor. To overcome these difficulties, comparative FE simulations are 
accomplished and calibrated to the experimental output. To obtain a common basis 
of comparison, the 3D results of the simulation have to be transferred into synthetic 
2D phase images. This post processing reproduces the physical transillumination 
according to the experimental photoelasticity. The combined numerical and exper- 
imental approach allows a reliable calibration of the FE model, providing verified 
3D stress fields as a basis for the presented studies. 

The paper is organized as follows: After a short description of the experimental 
procedure in section 2, section 3 addresses the FE model and simulation. Section 4 
covers the calibration and the interpretation of the experimental results. By means 
of the verified FE simulation, section 5 presents the composition of the stress tensor 
existing in the matrix before and after fibre fracture and its effect on the neighbour- 
ing fibre. Further studies focus on the effect of fibre-to-load orientation on both 
fibre and matrix loading and particularly on the shear component of the stress ten- 
sor along the fibre-matrix interface. Finally the progressive fracture behaviour at 
different load angles is investigated, followed by a discussion in section 6. 



2 Experimental photoelasticity with model samples 

Photoelasticity is a non-destructive method of observing preferably 2D stress dis- 
tribution in a specimen [9]. For that purpose a photoelastically active medium is 
required, such as polycarbonate Makrofol DE 1-1 (by Bayer [10]) used as ma- 
trix material. When transilluminating the specimen (see fig. 1), the linear polarized 
light is separated into two rectangular components according to the existing prin- 
cipal stress directions. Due to loading, the photoelastically active medium changes 
its refraction index, resulting in a phase shift between the two components of light, 
which have passed through the specimen at different velocities. Superposing both 
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components, erasement respectively duplication of amplitude is the consequence. 
The resulting projection of the specimen is a grey scale image including isochro- 
mates (lines of equal principal stress differences) and isoclines (lines of equal prin- 
cipal stress directions). Our experimental setup uses circular polarized light which 
is indifferent against direction. Therefore no isoclines appear. Further details on 
the experimental procedure can be found in [5,6]. The description of the follow- 
ing transformation into photoelastic phase images is published in [6]. The resulting 
images clearly indicate the darkest and lightest isochromates by a black-to-white 
phase jump containing information about the direction of the gradient. If the spec- 
imen's thickness d and the order of a phase shift n at a single point of the image 
is known, the existing principal stress difference o\ — a 2 can be calculated for the 
whole image using the photoelastic constant S according to 

<7i — (72 = S • n/d. (1) 

S characterizes the photoelastic activity. Since the laser beam passes the speci- 
men in reflection mode twice, d in Eq.l is twice the specimen thickness. Note 
that S behaves non-linearly for increasing load approximated by the experimen- 
tally determined function S = /(<7i — er 2 ) [9]. However, more accurately treated 
S — /(<7i, <7 2 ) is a function of o\ and a 2 . The final experimental output is a 2D 
field of integrated, scalar principle stress differences which means that neither in- 
formation on stress variations in specimen thickness nor components of stress are 
available. As the experiment deals with a 3D state of stress, the result does not 
reproduce the actual principal stress difference. It only accounts for those compo- 
nents existing in the plane rectangular to the way of light. For clarity, the result 
value is called secondary principle stress difference Aa sps . 

Both the experimental and numerical program contain pure matrix specimens to 
define the photoelastic constant S, single fibre specimens (fibre diameter ks 70 — 
80/xm) and specimens containing 5 parallel fibres oriented 0°, 5°, 10° and 15° to- 
wards the direction of load. Note that all stress indices throughout this work refer 
to fig. 1 with the triad aligned with the fibre orientation. 



3 The FE-simulation 

The model geometry corresponds to the used specimen in fibre length, diameter and 
distance, matrix thickness and positions of fibre fractures. The overall dimensions 
of the model were chosen smaller than those of the specimen, but large enough to 
observe far field stress. Fig. 1 includes the model dimensions in comparison to the 
specimen size [mm]. To save computation time, the model is divided into an outer, 
roughly discretized and an inner, finer discretized section embedding the fibres. 
The inner section (quarter-model) of a meshed five-fibre model (0°) including the 
cut-out for the embedded fibres is magnified in fig. 2. The irregularities of the mesh 
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Figure 1. Geometric situation of a five-fibre specimen (15°), including the two sections 
covered by simulation (dashed) and light path for reflection mode 

in the fibre plane are due to the demand, that the neighbouring matrix node has 
to be exactly in the plane of fracture to obtain comparable results at each fracture 
position. 




Figure 2. Quarter-model of inner section with 8-node thermally coupled brick elements 
with trilinear displacement and temperature using reduced integration 



Elastic-plastic material behaviour is used for the matrix description combined with 
linear-elastic fibres. Since the simulation covers the cooling-down process of the 
specimen during preparation, the material definition also contains linear thermal 
expansion coefficients for both materials. 
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4 Combined numerical and experimental approach 

4. 1 Synthetic phase images and calibration 

To obtain reliable FEM results, the material properties were calibrated using the 
experimental photoelastic phase images. To obtain a common basis of comparison, 
the FEM result is integrated in the direction of the model thickness. This post pro- 
cessing only accounts for the secondary principal stress difference Aa sps , which 
is the photoelastically relevant component of the stress tensor (see sec. 2). At each 
node Aa sps is obtained by a simple 2D principal axis transformation of the involved 
components of stress according to 

Aa sps = a 1 -a 2 = yj (a u - a 33 ) 2 + 4<x? 3 . (2) 

Due to the non-linear behaviour of the photoelastic activity (see sec.2), Aa sps has 
to be transformed into the corresponding order of phase shift n before integration. 
This value is calculated and normalized to unit thickness (d = 1) for each node 
using Eq. (1). Its weighted sums result in a 2D field of orders of phase shift n, 
which is visualized by a synthetic phase image [11]. These images are used for 
visual matching with the experimental phase images to determine the calibrated 
FEM material properties like the constitutive behaviour of the fibre material and its 
thermal expansion coefficient. Note that literature values for the PoiSSON ratio of 
the fibre and matrix material (y^ = 0.21, v m = 0.39) are used. 

The experimentally measured constitutive behaviour of the matrix was used for 
the FE simulations. Based on the implemented matrix behaviour, the fibre stiffness 
is calibrated matching the number of phase jumps of the experimental (fig. 3a) 
and synthetic (fig. 3b) phase images. Note that the relieve of the matrix by the 
embedded fibre leads to phase jumps surrounding the central part of the fibre. The 




Figure 3. Comparison of experimental (a) and synthetic (6) phase images of single-fibre 
specimen at a an applied tensile load of 54 MPa 
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linear thermal expansion coefficient of the matrix is calculated by the shrinkage 
of the specimen in direction of thickness during cooling-down. Furthermore, the 
resulting residual thermal stresses are used for the adaptation of the linear thermal 
expansion coefficient of the fibre. 



4.2 Interpreting the photoelastic results 



A direct advantage of the combined approach is the possibility to estimate the sig- 
nificance of the experimental photoelastic results. For example the FE result allows 
a comparison of experimental integral results with the actual local stress states. Fig. 
4 deals with a 0° single-fibre specimen after one fibre fracture located at 1.9mm 
with an applied tensile load of 27 MPa. The difference between the integral and 
the mid-plane value of Aa sps is an indirect measure for the homogeneity of the 
loading over specimen thickness and proves the necessity for 3D treatment. As this 
difference tends towards zero at the position of the fibre fracture, the increased ma- 
trix loading due to the fibre break extends over the whole specimen thickness. The 
large difference at less loaded positions is an evidence for a higher fluctuation of 
Aa sps along the light path. The definition of an algorithm to correct the deviations 




Figure 4. Comparison of Aa sps for a single-fibre specimen containing one fracture 



remains difficult because of the field of the photoelastically relevant Afx sps . In op- 
posite to the actual rotational symmetric state of stress around the fibre, Aa sps is 
not axisymmetric (see fig. 5). An additional challenge for formulating a correction 
algorithm is the necessary decomposition of the scalar photoelastic result. 
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Figure 5. Aa sps in plane of fibre break in a 0° single-fibre specimen at an applied tensile 
load of 56 MPa 



5 Results 



The presented FE results are based on the calibrated single fibre model with respect 
to material properties. This model is used for the fundamental characterization of 
the fibre-matrix interaction by showing the composition of the secondary principle 
stress difference (see sec. 5.1). In the following, the problem is extended to a five- 
fibre model for studies on load redistribution due to fibre fracture and interaction of 
neighbouring fibres (see sec. 5.2). Furthermore, results on the influence of fibre-to- 
load angle alteration on the fibre loading (see sec. 5. 3.1), the loading of matrix (see 
sec.5.3.2) and the progressive fracture behaviour (see sec.5.3.3) are presented. 



5. 1 Composition of the secondary principal stress difference 



Due to the calculated stress tensor the contribution of the different stress compo- 
nents to the photoelastically visualized Aa sps is known (see fig. 6). It is clearly 
visible that the main component is the stress in fibre direction a 33 . The absolute 
value of the shear stress <r 13 has its maxima at both ends of the fibre, due to the 
load transmission from matrix to fibre. The cooling-down process of the specimen 
during preparation leads to a compressive stress component an around the fibre. 
Due to the rather complex composition of A<r sps even for single-fibre specimens 
aligned with the load direction, the interpretation of the scalar experimental output 
seems to be difficult. 
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Figure 6. Components of Aa sps in the matrix along fibre-matrix interface in a 0° sin- 
gle-fibre specimen without fracture at an applied tensile load of 25 MPa 

5.2 Load redistribution due to fibre fracture 



Load redistribution due to a fibre fracture leads to an increased matrix loading in 
a butterfly-shaped area around the fracture (see fig. 12a and earlier photoelastic 
investigations [5]). This influences both, the state of stress of the matrix along the 
neighbouring fibre-matrix interface and the fracture behaviour of the neighbouring 
fibre itself. It is reflected in a remarkably high probability of the fibre to break in the 
sphere of influence of a neighbouring fracture. Two pairs of consecutive fractures 
are shown in fig. 7. The FE simulation proves the butterfly shape of Aa sps . 




Figure 7. Synthetic phase image of 5-fibre 0° specimen after 5 fibre fractures at an applied 
tensile load of 45 MPa 

For a quantitative description of the load redistribution, Aa sps and its components 
are extracted along fibre-matrix interfaces (see fig.8). Note that fibre labels refer 
to fig. 7. Looking at the direct surrounding of the breaking of fibre 4, the fracture 
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Figure 8. a: Aa sps in matrix before and after fracture of fibre 4; b: Components of Aa sps 
at the neighbouring fibre-matrix interface 3 before and after fibre fracture (thin/bold lines) 

causes a maximum of Aa sps in the matrix (see fig. 8a) consisting mainly of two 
stress components: stress in fibre direction cr 33 and shear stress <7i 3 . The increase 
of 0-33 is caused by the complete tensile load redistribution from fibre to matrix, 
while <7i3 increases because of two newly emerged areas of load transmission from 
matrix to fibre. 

The consequences on the neighbouring fibre 3 are shown in fig. 8b where the frac- 
ture leads to a minimum of Aa sps in the fracture plane with two surrounding max- 
ima. This typical shape is a cut through the butterfly area. Furthermore, fig. 8a 
includes the decreasing influence of the existing fracture in fibre 2 on the loading 
along fibre 3 and 4 (i) and the influence of the fracture of fibre 5 on fibre 4 (ii) be- 
fore it consecutively breaks at the increased local maximum of Aa sps (iii). Looking 
at the composition of Aa sps in an enlarged zone along the neighbouring fibre 3 (see 
fig. 8b), (T33 increases in the fracture plane analogously to the direct fracture sur- 
roundings in fibre 4. At the same position, the absolute value of compressive stress 
(Tn decreases to about MPa, which might weaken the fibre-matrix bonding. The 
maximum of the absolute value of a 13 causes the two local maxima of Aa sps , which 
agree with the experimentally observed cascade positions. The butterfly shape and 
the probability of the neighbouring fibre to break in typical cascade formation (see 
fig. 7, fibres 5 and 4) therefore arise from the behaviour of the shear stress <ti 3 . The 
comparable relevance of 033 for the damage evolution is addressed in sec.5.3.3. 

5.3 Influence of fibre orientation 

To study the progressive failure of long fibre composites under shear stress loading, 
the test series contains several 5-fibre specimens at an increasing load angle. By 
means of the presented simulations we can extract the states of stress during the 
experiment, giving rise to interpretation of the gradual failure. 
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5.3.1 Influence of fibre orientation on the loading of fibre 



Since photoelastic experiments show the effect of the embedded fibre but not the 
stress state of the fibre itself, the investigation relies on the results of FE simula- 
tions. The presented study focuses on the tensile stress <7 33 , which is assumed to 
be one reason for fibre failure. 0-33 increases rapidly in the areas of load transmis- 
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Figure 9. Tensile stresses in the outer fibres of the five-fibre package at an applied load of 
30 MPa 

sion up to a relatively constant value over the mid-fibre positions. Looking at the 
influence of variation of fibre orientation, the expected decrease of load absorption 
by the fibres with increasing load angle is observed and in agreement with experi- 
ments. Additionally the numerical results prove the experimental observation, that 
outer fibre tends to break first. The tensile stress in the outer fibre of a package is 
about 5% higher than in the middle fibre. 



5.3.2 Influence of fibre orientation on the loading of matrix 

To understand the influence of load angle alteration on the stress distribution exist- 
ing in the matrix, Aa sps is decomposed with special interest in a 13 . Therefore the 
stresses along the fibre-matrix interface are chosen for further evaluation (see figs. 
10, 1 1). At higher load angles, the decreasing load absorption by the fibres results 
in a moderate increase of Aa sps (see fig. 1 1). Both the synthetic phase images (see 
fig. 10) and the extracted numerical output (see fig. 1 1) show a maximum of Aa sps 
in the area of load transmission at the right end of the fibre. Additionally, the max- 
imal relief of the matrix by the embedded fibre is shifted to the left. Looking at 
the components of Aa sps the stress in fibre direction 033 decreases slightly due to 
the rotation of coordinates but shows little asymmetric behaviour. In return an is 
shifted to a tensile stress (see for 15° fig. 13a]. The largest sensitivity against the 
load angle is found for the shear stress a 13 . Loaded at an angle, the shear stress is 
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Figure 10. Synthetic phase image of a five-fibre model at 15° load angle 




Position along fiber LtnmJ 



Figure 11. Components of Aa sps along the path of fig. 10 at 34 MPa 

a superposition of two components: The first one arises from the above-described 
process of load transmission from the matrix to the fibre (see sec. 5.1) leading to the 
extreme values at the fibre ends. The second component is constant, which origi- 
nates from the rotation of the fibre-matrix interface and shifts the shear stress cr 13 up 
to higher values. As a result, the inflection point in mid-fibre position is no longer 
at zero but describes the amount of the second component of cr 13 . 



5.3.3 Influence of fibre orientation on the fracture behaviour 

In order to understand the progressive failure of multi-fibre specimen, one needs 
to understand the influence of fibre orientation on the loading of fibre and matrix. 
First we focus on the direct fracture zone (see fig. 12b*), where those areas of the 
fibre-matrix interface, which form an acute angle with the fracture plane are higher 
loaded (see figs.l2b m ,13a m ). The similar effect was previously observed for the 
fibre ends (see sec. 5. 3. 2). Along the opposite interface a distinctive minimum of 
Aa sps is observed (see figs. 12b", 13a w ). The consequences of different fibre ori- 
entations can be seen in the numerically determined stress state (see fig. 13a). As 
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Figure 12. Photoelastic phase images of fracture areas. Thin lines define the fracture plane, 
while bold lines connect the points of largest Aa sps at increasing distance from the fracture 

shown in fig. 1 1 , <t 13 is again affected by the variation of fibre orientation and forms 
the crucial factor for the described minimum of Aa sps (see figs.l2b M , 13a M ). The 
numeric results also show that a u and <7 33 are almost insensitive with respect to the 
load angle, as followed from their symmetry to the fracture plane (see fig. 13a), and 
behave similarly to 0° specimens. Note that o" n is simply shifted to positive values 
due to a superposed tensile stress component exceeding the compressive residual 
stress around the fibre (see fig. 13a). 




Figure 13. a: Components of Aa sps after fibre fracture along the lowest fibre-matrix inter- 
face of fig. 12b at an applied load of 44 MPa and 15° b: Loading of the neighbouring fibre 
and the matrix along its interface 

Analyzing the surrounding matrix area, we observe that the lines of maximum 
Aa sps do not rotate according to the fibre orientation (compare figs. 12a and 12b), 
so the orientation of the 'butterfly wings' is not affected by the fibre orientation. In 
contrast, the value of Aa sps along the bold line is subject to the variation of fibre 
orientation, in such a way that the 'butterfly wings' are distorted asymmetrically. 
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Numerous experimental phase images including earlier investigations [6,8] reveal 
two different types of progressive fracture events: cascade-like propagation of fibre 
fractures (see fig. 7 fibres 4/5) and those in fracture plane (see fig. 7 fibres 1/2). The 
differences in the stress composition at the two predetermined fracture positions 
(fig. 13b:i, ii) are used for interpretation. The cascade-like propagation (see fig. 13b l ) 
is not due to the elevation of the tensile stress of the fibre a[ 3 but might be the result 
of the influence of the shear stress component cx^ at the fibre-matrix interface. On 
the other hand the maximum of a( 3 is the crucial factor for the second type of 
progressive fracture events and coincides with the maximum of a r 33 in plane of 
fracture (see fig.l3b M ). 

Including the influence of load angle variation, the neighbouring fibre has jsut one 
predestined cascade position of break, whereas at 0° the undistorted butterfly has 
two equally probable fracture positions per neighbouring fibre. Independent of the 
load angle, the fracture position is located on the bold line (see fig. 12b). 



6 Discussion 



As a basis for this work we have documented and evaluated the influence of the 
load direction on the progressive failure of multi-fibre specimens via confocal laser 
scanning photoelasticity (see fig. 12). The following calibration of the FE model on 
simple specimens containing one fibre gives the necessary material properties for 
the more sophisticated simulations (see sec.4.1). Extensive agreement of experi- 
mental and synthetic phase images allows reliable quantitative statements on the 
load redistribution during the failure process. In addition to the experimental pho- 
toelasticity, the FEM simulation provides a stress tensor used for further studies 
(see sec. 5.1). In principle the composition of the photoelastically visualized sec- 
ondary principle stress differences strongly depends on the applied tensile load and 
load orientation (see fig .11). We have shown that the shear stress component is an 
important factor e.g. for the butterfly shape of increased matrix loading due to fibre 
fracture. The significance of the shear stress increases for higher load angles (see 
sec.5.3.2). 

The FEM simulation includes the photoelastically inactive fibres and proves the 
expected reduction of load absorption by the fibres for increasing load angle (see 
sec. 5. 3.1). Earlier photoelastic investigations and the presented simulation results 
reveal two different types of progressive fracture events, one influenced by the shear 
stress component (cascade-like failure) and the other one by the tension in fibre di- 
rection (see sec. 5. 3. 3). Since the butterfly is distorted for angular load, only one 
position at the neighbouring fibre becomes the predestined cascade-like fracture 
position. Note that the orientation of the cascade is independent of the fibre direc- 
tion. 
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Despite of the simplifications and assumptions made in the model, the presented 
simulations contribute to the comprehension of load redistribution processes, the 
load angle influence, the fibre interaction and the fracture propagation. Particularly 
due to the iterative calibration process, it was necessary to neglect creep rates and 
adaptive fibre-matrix debonding. The presented results show, that the interpretation 
of photoelastic phase images visualizing complex states of stress has to be broad- 
ened by a combined numerical investigation to gain full insight. Furthermore, the 
presented numerical approach can be extended beyond the physical limits of photo 
elasticity such as small fiber distances or complex 3D fibre arrangements. 
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